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Dr.  George  Yoder,  Code  1131 
Office  of  Naval  Research 
Materials  Division 
800  North  Quincy  Street 
Arlington,  VA  22217 

Dear  George: 

Th'S  letter  is  meant  to  be  a  progress  report  for  the  period  March  1, 1992  to  February  28, 1993 
of  our  grant  with  ONR  #N000l4-92-J-1752.  I  am  also  enclosing  a  PPPH  report  for  your  records. 

Although  this  is  the  first  year  of  this  proposal,  we  have  already  made  significant  progress 
towards  our  goal  to  be  able  to  accurately  theoretically  model  and  understand  the  strengths  of 
Matrix/Second-phase-particle  Interfaces.  Our  accomplishments  have  been  along  two  directions. 
The  first  involves  using  our  state-of-the-art  ab-initio  dynamics  approach  to  calculate  the  properties 
of  the  bare  Cu(100)  surface.  This  is  needed  as  a  test  before  we  may  begin  to  study  the  interfaces 
with  the  second-phase-particles.  The  results  of  our  calculations  of  some  structural  properties  are 
summarized  in  Table  I. 

Table  I 

Comparison  between  theory  and  experiment  of  various  ground  state 
properties  of  bulk  Cu  and  its  (100)  surface 


Properties  of  Cu 

Theory 

Experiment 

Lattice  constant  in  bulk  (A) 

3.63 

3.60 

Bulk  modules  (M  Bu) 

1.56 

1.42 

Bulk  cohesive  energy  (eV/atom) 

3.91 

3.50 

Surface  layer  relaxation 

-4% 

-3% 

Surface  energy  (eV/atom) 

0.84 

0.74  ±  .2 

As  you  can  see  the  agreement  between  theory  and  experiment  is  excellent.  This  is  particularly 
true  when  you  consider  that  no  experimental  input  enters  the  theory  except  for  the  atomic  number 
of  Cu.  The  electronic  charge  density  of  the  (100)  surface  is  shown  at  the  top  panel  of  Fig.  1.  In 
order  to  get  a  better  understanding  of  the  effects  of  the  s-electron  and  d-eleciron  contributions  to 
this  density  we  show  at  the  bottom  of  Fig.l  the  density  obtained  by  subtracting  off  the  atomic  d- 
electrons.  The  s-electrons  are  seen  to  be  very  nearly  uniformly  spread  out,  almost  as  if  they  were 
free-electron  like.  This  feature  could  be  used  in  modeling  the  region  away  from  the  Cu  interface  as 
a  jellium  system. 

As  far  as  the  details  of  the  interface  are  concerned,  we  have  discovered  that  we  can 
theoretically  design  a  Cu(100)  -  TiC>2  interface  that  is  lattice-matched  to  less  than  0.5% !  The 
interface  supercell  is  Vl3  x  Vl3  for  Cu  and  2x2  for  Ti02-  This  system  will  than  provide  a 
tractable  and  interesting  interface  to  study  using  our  ab-initio  techniques,  and  is  also  possibly  an 
interface  that  could  be  made  experimentally.  We  are  consistently  in  contact  with  Professor  Warren 
Garrison  and  his  group  at  Camegie-Mellon  and  have  suggested  this  new  interface  as  an  interesting 
system  to  study.  Professor  Garrison  has  recently  finished  some  interesting  studies  of  secondary 
crystallization  and  most  recently  has  gone  on  to  prepare  materials  that  are  internally  oxidized.  It 
appears  that  both  theory  and  experiment  are  now  at  the  point  to  begin  studies  of  the  interphase 
interfaces. 

The  second  accomplishment  this  past  year  involves  the  introduction  of  a  new  theoretical 
method  that  appears  to  be  extremely  promising  for  performing  ab-initio  calculations  of  complex 
material  systems.  Although  the  technique  is  still  in  a  very  preliminary  form,  and  not  yet  capable  of 
studies  of  the  types  of  material  systems  of  interest  for  this  proposal,  it  has  great  potential.  The 
technique  involves  the  introduction  of  Wavelet  Theory  and  Multi-resolution  Analysis  to  the 
electronic  structure  ajena.  Test  calculations  on  reproducing  the  Is  states  of  all  the  atoms  in  the 
periodic  table,  as  well  as  the  ground  state  of  the  hydrogen  molecule  ion,  successfully  demonstrate 
the  potential  power  of  this  new  technique.  The  results  for  the  hydrogen  molecule  ion  are  show  in 
Fig.  2.  Here  the  open  circles  represent  tne  exact  solution.  We  have  already  written  a 


paper  on  this  subject  and  sent  it  out  for  publication,  a  copy  of  which  is  enclosed  with  this  letter. 
We  are  very  excited  about  the  prospects  of  this  new  approach,  and  plan  to  pursue  it  vehemently. 


I  hope  this  is  sufficient.  If  you  need  any  additional  information  please  let  me  know. 


Sincerely, 


John  D.  Joannopoulos 
Professor  of  Physics 


JDJ:ia 


^r  ao  /  1^3 


Wavelets  in  Electronic  Structure  Calculations 

K.  Cho,  T.  A.  Arias,  and  J.  D.  Joannopoulos 
Department  of  Physics,  Massachusetts  Institute  of  Technology,  Cambridge,  MA  02139 

Pui  K.  Lam 
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A  three  dimensional  (3D)  wavelet  analysis  is  employed  to  develop  a  new  formal¬ 
ism  for  electronic  structure  calculations.  The  wavelet  formalism  provides  a  systemati¬ 
cally  improvable  and  tractable  description  of  electronic  wavefunctions  and  overcomes 
limitations  of  conventional  basis  expansions.  The  potential  power  of  the  wavelet 
formalism  for  ah  initio  electronic  structure  calculations  is  demonstrated  by  a  calcu¬ 
lation  of  Is  states  for  all  the  naturally  occurring  nuclei  on  the  periodic  table  and  the 
interaction  energies  of  the  hydrogen  molecule  ion. 
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Spectral  analysis  using  the  Fourier  transform  (FT)  is  a  powerful  method  for  solv¬ 
ing  many  problems  in  science  and  engineering.  This  method,  however,  is  not  appro¬ 
priate  for  problems  that  require  a  localized  description  in  real  space  and  in  Fourier 
space.  Although  several  techniques,  such  as  the  windowed  Fourier  transform,  have 
been  invented  to  attempt  to  overcome  this  limitation,  by  far  the  most  promising  new 
technique  is  the  wavelet  transform  (WT)  [1].  Unlike  the  plane  wave  (PW)  basis  func¬ 
tions  of  the  FT,  the  wavelet  basis  functions  are  localized  both  in  real  space  and  in 
Fourier  space  [1-4].  Furthermore,  multiresolution  analysis  (MRA)  of  the  WT  pro¬ 
vides  a  systematic  successive  approximation  scheme  for  practical  applications  [1,2]. 
Applications  of  the  WT  have  focussed  primarily  on  digital  signal  processing  (ID), 
compact  image  coding  (2D),  and  related  fast  mathematical  algorithms  [3 j .  Most  re 
cently  the  WT  has  been  applied  to  the  analysis  of  chaotic  behavior  and  turbulence 
in  2D,  the  coherent  states  of  quantum  optics  and  quantum  field  theory,  and  to  real 
space  renormalization  group  theory  [4]. 

In  this  letter,  we  introduce  a  new  method  for  the  application  of  WT,  in  31).  to 
electronic  structure  investigations  of  material  properties.  The  traditional  ab  initio 
total  energy  density  functional  pseudopotential  methods  for  electronic  structure  cal¬ 
culations  use  either  LC AO-type  (c.g.  atomic,  Gaussian,  or  Slater)  basis  sets  or  PW 
basis  sets  [5 j .  LCAO-type  basis  sets  typically  are  capable  of  describing  the  electronic 
structure  with  a  small  number  of  basis  functions,  but  lack  an  explicit  and  well-defined 
procedure  for  systematic  improvement  [6],  and  generally  complicate  the  calculation  of 
forces  because  of  the  presence  of  Pulay  terms  ] 7] .  On  the  other  hand,  PW  basis  sets 
provide  a  systematic  expansion  of  electronic  wavcfunctions,  and  simple  algorithms  for 
calculating  force*,  but  sufler  from  the  same  limitations  as  the  FT  me! hod.  Specifi 
rally,  the  description  of  the  rapid  variations  of  the  electronic  wavefunctmn  close  to 
the  atomic  nuclei  (atomic  core)  requires  a  large  number  of  PW  functions  even  though 


the  volume  of  the  atomic  core  is  only  a  small  part  of  the  total  volume  of  the  sys¬ 
tem.  In  addition  PW  basis  sets  are,  in  prin  iple,  not  well  suited  for  the  description 
of  isolated  molecules  and  surfaces.  This  is  because  the  periodic  nature  of  the  basis 
introduces  unnecessary  periodic  images  and  a  redundant  high  resolution  description 
of  the  vacuum  regions. 

The  “best  of  all  worlds”  then  would  be  to  have  a  basis  set  that,  unlike  a  PW 
basis,  can  zoom  into  the  atomic  core  regions  and  does  not  introduce  unnecessary 
periodic  images,  and,  unlike  an  LCAO-type  basis,  does  not  have  Pulay  terms  and 
provides  an  explicit  prescription  for  complete  expansion  of  the  electronic  wavefunc- 
tions.  The  wavelet  basis  has  all  these  desired  properties:  multiscale  decomposition  ol 
the  wavefunction,  the  localized  description  of  the  rapid  variations,  and  systematically 
complete  expansion  similar  to  the  PW  expansion.  Therefore,  as  we  demonstrate  in 
this  letter,  the  wavelet  basis  overcomes  the  limitations  of  both  conventional  basis  sets 
and  retains  only  their  advantages. 

We  begin  our  discussion  with  a  brief  review  of  the  basic  concepts  underlying 
wavelets.  We  then  extend  the  wavelet  formalism  to  the  solution  of  the  Schrbdinger 
equation  with  the  introduction  of  spherically  symmetric  basis  functions.  As  a  simple 
example  we  apply  our  technique  to  the  study  of  a  hydrogen  atom.  We  then  demon¬ 
strate  that  a  single  small  basis  set  is  capable  of  calculating  accurately  the  Is  states 
of  all  the  nuclei  from  Hydrogen  through  Uranium!  Finally  we  demonstrate  that  a 
wavelet  basis  can  easily  describe  covalent  bonding  and  illustr.  *e  its  use  with  the  hy¬ 
drogen  molecule  ion  (Hf  ).  We  conclude  with  a  discussion  of  the  straightforward 
extension  of  the  current  analysis  to  periodic  systems  and  all  electron  calculations. 

The  wavelet  transform  and  multiresolutwn  analysis.  -  Given  a  square  integrable 
function  space  L2(R 3),  wavelets  impose  a  hierarchical  structure  of  subspaces  with 
different  resolutions  which  forms  a  mult irrsolut  ion  analysis.  1'he  space  of  functions  at 


resolution  2~J  (or  more  simply  at  resolution  j)  is  represented  by  l'\ ( R3 ) ,  and  spanned 
by  the  basis  set  of  scaling  functions  at  resolution  j,  {#j.n(r)},  where  n  specifies  the 
center  of  the  basis  function.  The  hierarchical  structure  is  then 

-  V.2  C  C  l'o  C  Vj  C.  Vj  ■■■■  fl) 

I  he  approximation  space  Vj  is  decomposed  into  a  sum  of  the  coarser  approximation 
space  Vj_i  and  the  wavelet  space  W which  describes  the  detail  at  resolution  j 
(Vj  —  Vj_ J  ©  WO-,)  and  is  spanned  by  the  basis  set  of  the  wavelet  functions  at 
resolution  j,  {V’j,n(r)}-  With  a  repetition  of  this  decomposition,  the  L7{  R 3)  space 
can  be  expressed  as  either  the  sum  of  the  wavelet  spaces  of  all  resolutions  or  the  sum 
of  one  approximation  space  and  the  wavelet  spaces  of  higher  resolutions: 

0  Wj  -  Vjo  ©  0  Wj  =  r,2(R3).  (2, 

J  }=jo 

Therefore,  any  square  integrable  function  /( r)  can  be  expanded  either  as  a  sum  of  the 
wavelet  functions  of  all  resolutions  or  as  a  sum  of  the  scaling  functions  at  resolution 
j  —  jo  and  the  wavelet  functions  of  all  finer  resolutions  j  >  j0-  Iti  this  work,  we 
will  use  the  latter  expansion  because  it  introduces  approximations  only  at  the  high 
resolution  cutoff  in  practical  applications. 

With  the  introduction  of  two  projection  operators,  A3  (approximation  operator) 
and  Dj  (detail  operator),  which  project  a  function  into  Vj  and  11  j  respectively,  one 
may  express  /( r)  as 

/(r)  -  AJ0f(r)  +  Y,  bjf(r).  (ill 

J  =  JO 

The  approximation  and  the  detail  of  a  function  /( r)  at  resolution  ;  are  expanded  in 
terms  of  the  basis  functions, 

Ajf(r)  —  ^  ^  ^i.u^uijr).  !  I ) 
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Combination  of  Eqns.  (3)-(5)  and  truncation  at  the  finest  resolution  }mar  leads  to 
the  wavelet  expansion  of  /( r)  as 

__  Jmar 

/(r)  =  Eajo.n^.,o.n(r)  +  n(r).  (6) 

«  3  =  3  0  n 

This  expansion  still  contains  an  infinite  number  of  basis  functions  associated  with 
the  lattice  {n}.  Since  the  scaling  functions  and  the  wavelet  functions  are  spatially 
localized,  one  may  retain  only  the  basis  functions  that  have  significantly  large  coeffi¬ 
cients  in  Eqn  (6)  for  the  problem  at  hand.  This  leads  to  an  expansion  with  a  finite 
number  of  basis  functions  and  allows  one  to  use  different  resolutions  for  different  spa 
tially  localized  regions.  In  particular,  for  the  description  of  the  rapid  oscillations  of 
the  electronic  wavefunction  in  an  atomic  core  region,  one  need  add  higher  resolution 
scales  only  to  the  core  region  and  thereby  systematiclly  improve  the  calculation. 

Construction  of  a  wavelet  basis  set.  -  Although  the  wavelet  and  scaling  functions 
are  far  from  unique,  we  found  that  the  following  forms  are  particularly  convenient  in 
practical  applications.  The  wavelet  functions  r)  are  chosen  to  be  the  Mexican 
hat  functions  (the  Laplacians  of  Gaussians)  which  form  a  fairly  tight  frame  !lj.  The 
scaling  function  4>o.n(  ^corresponding  to  the  Mexican-hat  function  is  chosen  to  be  a 
Gaussian  function.  Both  the  Mexican-hat  functions  and  the  Gaussian  functions  are 
spherically  symmetric.  The  relationship  between  the  Gaussian  scaling  function  and 
the  Mexican-hat  wavelet  function  is  not  exact,  but  is  quite  a  good  approximation 
as  discussed  in  Ref.  1.  These  localized  basis  functions  are  spatially  arranged  so 
that  their  centers  form  a  simple  cubic  lattice  [8j.  This  is  shown  schematically  in 
the  top  panel  of  Fig.  1.  The  lattice  spacing  c/0  at  resolution  j  -  0  is  chosen  small 
enough  to  give  a  fairly  tight  frame  j  1 ,9 j .  The  basis  functions  with  resolution  j  are 
arranged  correspondingly  on  the  lattice  sites  with  lattice  spacing  </0/ 2-’.  Ihe  centers 
of  basis  functions  for  different  resolutions  are  selected  so  as  not  in  overlap,  and  the 


centers  of  all  the  basis  functions  form  a  simple  cubic  lattice  with  spacing  do/ 2Jm“' .  In 
order  to  construct  a  finite  basis  set  we  collect  the  most  important  basis  function*  l>y 
introducing  spheres  of  finite  support  radii  at  each  nuclear  center.  This  is  illustrated 
schematically  at  the  bottom  panel  of  Fig.  1.  Correspondingly  smaller  support  radii 
are  chosen  for  correspondingly  higher  resolutions  so  that  deep  core  regions  have  more 
resolution  scales.  We  designate  this  finite  basis  set  as  {<£0.n(r),  V'J,n{r)'<  Jma*  >  )  >  0} 
or  more  compactly  as  {f>r(r)}  where  i  is  a  simplified  notation  for  {j,n}. 

Solution  of  Schrodinger  equation.  -  The  Schrodinger  equation,  >  -  cf'k  >, 

for  electronic  structure  calculations  becomes  a  generalized  eigenvalue  problem  in  a 
general  nonorthogonal  wavelet  basis.  Expansion  of  the  wavefunclion  in  terms  of  the 
basis  function  (I'k  >=  c^fe,  >)  leads  to  the  secular  equation, 

E  <  b,\H |6,  >  c,  =  <  bfib,  ;■  c„  (7) 

*  I 

With  our  choice  of  wavelet  and  scaling  functions,  .11  the  matrix  elements  in  Kq.  (7) 
can  be  calculated  analytically.  Solution  of  Eq.  (7/  may  then  proceed  by  a  number 
techniques  including  square  root  matrix  diagonabza.ion  [10],  Oholesky  decomposi- 
tion  [11],  molecular  dynamics  [5],  and  conjugate  gradients  |5j  approaches.  In  this 
work  we  have  used  both  square  root  matrix  diagonalization  and  conjugate  gradients 
approaches. 

The  hydrogen  atom.  -  In  order  to  gain  a  sense  of  the  optimal  values  of  the  various 
parameters  of  the  formalism,  we  chose  to  study  the  hydrogen  atom  as  a  simple  test 
case.  After  performing  calculations  using  many  different  resolution  scales  and  different 
support  radii,  we  find  ^at  a  satisfactory  minimal  basis  set  for  the  Hydrogen  atom 
contains  25  basis  functions  with  four  resolution  scales:  (7  scaling  function*  lor  /  0. 

plus  6  wavelets  for  j  =  0,  plus  6  wavelets  ior  j  -  1 ,  plus  (i  wavelets  fur  /  2)  which 

we  designate  simply  as  ( 7  |  three  fi’s).  These  correspond  respectively  to  support 
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radii  of  lau,  0.5au,  0.2bau,  and  0.125au.  The  calculated  Is  eigenenergy  is  within  2% 
of  the  exact  value,  and  the  calculated  Is  radial  wavefunction  is  shown  in  Fig.  2  (top 
panel).  The  calculated  wavefunction  ditiers  from  the  exact  one  primarily  near  the 
origin  (within  the  resolution  limit  jmax)-  The  sinali  difference  at  r  —  2au  is  due  to 
the  finite  support  radii  used  for  the  basis  set.  By  increasing  the  support  radii,  the 
difference  at  r  =  2 au  disappears  as  shown  in  the  bottom  panel  of  Fig.  2.  The  result 
in  the  bottom  panel  of  Fig.  2  corresponds  to  a  basis  set  of  85  basis  functions  with 
three  resolution  scales  (33  -f  two  26’s)  with  support  radii  (2gu,  lau,  and  0.5au),  and 
gives  the  Is  eigenenergy  to  within  0.5%  of  the  exact  value. 

Hydrogen  to  Uranium.  -  By  adding  more  resolution  scales,  one  can  calculate 
the  wavefunctions  of  heavier  nuclei  within  the  same  accuracy  as  follows.  As  the 
atomic  number  increases,  the  Coulomb  potential  becomes  stronger,  and  consequently 
higher  resolution  scales  are  needed  in  the  core  region.  However,  only  one  additional 
resolution  scale  needs  to  be  added  each  time  the  atomic  number  doubles.  Therefore, 
for  the  description  of  Is  wavefunctions  from  ffydrogen  to  Uranium,  one  needs  <o  use  a 
basis  set  with  eleven  resolution  scales  (7+  ten  6’s).  Using  this  fixed  basis  set  consisting 
of  only  67  basis  functions,  one  can  calculate  (to  within  3%)  the  Is  eigenvalues  of  all 
the  naturally  occurring  nuclei  on  the  periodic  table!  The  results  are  shown  in  Fig.  ,5 
and  include,  for  simplicity,  only  the  even  nuclei.  Note  that  the  eigenvalue  for  each 
nucleus  is  larger  than  the  exact  vaiue  by  very  nearly  a  constant  percentage  so  that 
the  results  appear  to  lie  on  a  straight  line  as  expected. 

The  hydrogen  molecule  ion.  ■  We  new  investigate  the  efficiency  of  the  wavelet  basis 
for  the  description  of  chemical  bonds.  As  a  simple  example,  we  consider  the  energy 
of  an  //2+  molecule  as  a  function  of  the  separation  H  between  the  two  protons.  The 
total  energy  (electronic  eigenenergy  plus  the  Coulomb  potential  energy)  is  plotted  as 
a  function  of  R  in  Fig.  -1.  The  basis  set  for  the  calculation  is  selected  using,  the  same 
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support  radii  as  for  the  hydrogen  atom  calculation  with  the  larger  basis  set  (Fig.  2 
bottom  panel),  and  the  total  number  of  basis  functions  varies  between  141  and  167 
[12]  depending  on  R.  It  is  important  to  note  that  the  centers  of  the  basis  functions 
are  not  following  the  nuclei  as  the  separation  R  changes;  i.c.,  the  underlying  infinite 
basis  set  is  not  correlated  to  the  ionic  positions.  Rather,  it  is  the  support  spheres, 
which  are  attached  to  each  nucleus,  that  designate  and  determine  the  finite  basis  set 
as  illustrated  in  the  bottom  panel  of  Fig.  1.  This  property  of  the  basis  set  does  not 
lead  to  Pulay  forces,  and  allow's  the  wavelet  functions  to  overcome  the  limitations  of 
LCAO-type  basis  functions  [7].  Finally  we  note  that  the  calculated  values  in  Fig.  4 
are  larger  than  the  exact  values  (13)  by  a  nearly  constant  amount  of  1%.  This  implies 
that  the  wavelet  basis  also  gives  an  excellent  representation  of  ionic  forces. 

Conclusions.  -  We  have  developed  a  3D  wavelet  formalism  that  is  applicable  in 
general  to  electronic  structure  investigations  of  materials,  and  have  demonstrated  its 
potential  power  by  performing  calculations  of  the  deepest  core  stales  of  all  the  nuclei 
on  the  periodic  table  and  the  interaction  energy  of  an  molecule.  The  wavelet 
expansion  with  a  small  number  of  basis  functions  gives  an  accurate  description  of 
the  wavefunction  both  in  an  atomic  core  region  and  in  a  chemical  bond.  To  extend 
this  approach  to  many-electron  systems  within  LDA,  one  can  simply  use  traditional 
techniques  for  inclusion  of  Hartree  and  exchange-correlation  interactions.  In  partic¬ 
ular,  schemes  developed  for  localized  basis  sets  [14]  provide  one  natural  framework 
for  this  extension.  Finally,  the  extension  of  the  current  analysis  to  periodic  systems 
is  straightforwardly  achieved  by  introducing  a  Bloch  transformation  of  the  wavelet 
functions  [15].  With  the  extensions  to  many-electron  systems  and  periodic  systems, 
the  wavelet  formalism  can  open  a  completely  new  direction  of  the  development  for  ah 
initio  total  energy  calculations. 
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FIGURES 


FIG.  1.  Top  panel  shows  a  two  dimensional  arrangement  of  the  centers  of  basis  func¬ 
tions  on  a  simple  square  lattice.  Open  circles  are  the  centers  of  the  basis  functions  at 
resolution  j  —  0,  and  filled  circles  are  the  centers  of  the  basis  functions  at  resolution  j  -  1. 
Bottom  panel  shows  the  basis  functions  that  lie  within  each  support  radius  that  is  centered 
on  an  atomic  position  (diamond-symbol).  The  larger  (smaller)  radius  corresponds  to  j  —  0 
(j  =  1)  resolution. 


FIG.  2.  Top  panel  shows  the  Is  radial  wavefunction  of  a  hydrogen  atom  calculated 
with  25  basis  functions  (7  +  three  6’s).  Bottom  panel  shows  the  Is  radial  wavefunction  of 
a  hydrogen  atom  calculated  with  85  basis  functions  (33  +  two  26’s).  The  continuous  lines 
are  from  the  wavelet  calculations,  and  the  broken  lines  are  the  exact  wavefunctions. 


FIG.  3.  The  filled  circles  are  the  Is  eigenenergies  of  all  the  nuclei  with  even  atomic 
numbers  on  the  periodic  table  (Z  =  2, 4,. ..,92)  calculated  with  the  fixed  67  basis  function 
set. 


FIG.  4.  The  total  energy  of  a  molecular  hydrogen  ion  is  shown  as  a  function  of  the 
separation  between  the  nuclei.  The  filled  circles  are  the  calculations  with  141  to  167  basis 
functions  depending  on  the  separation  Ft.  The  open  circles  are  the  exact  values. 
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